library(ggplot2)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(tidyr)
library(gridExtra)
##
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
##
## combine
library(Hmisc) #for correlation matrix
## Loading required package: lattice
## Loading required package: survival
## Loading required package: Formula
##
## Attaching package: 'Hmisc'
## The following objects are masked from 'package:dplyr':
##
## src, summarize
## The following objects are masked from 'package:base':
##
## format.pval, units
NA.month <- read.csv("ClimateNA/outputs/monthly_tidy.csv")
inst.month.P <- read.csv("Instrumental/outputs/monthly_precipitation.csv")
inst.month.T <- read.csv("Instrumental/outputs/monthly_temperature.csv")
NA.month$Month <- as.factor(NA.month$Month)
inst.month.P$month <- as.factor(inst.month.P$month)
inst.month.T$month <- as.factor(inst.month.T$month)
m.m.t <- NA.month %>% select(Year, Month, Tave)
i.m.t <- inst.month.T %>%
select(year, month, mean_temp_year) %>%
filter(year > 1900) %>%
rename(Year = year,
Month = month,
Tave = mean_temp_year)
m.m.t$Source <- "ClimateNA"
i.m.t$Source <- "Instrumental"
months.T <- rbind(m.m.t, i.m.t)
m.m.p <- NA.month %>% select(Year, Month, PPT)
i.m.p <- inst.month.P %>%
select(year, month, total_precip_year) %>%
filter(year > 1900) %>%
rename(Year = year,
Month = month,
PPT = total_precip_year)
m.m.p$Source <- "ClimateNA"
i.m.p$Source <- "Instrumental"
months.P <- rbind(m.m.p, i.m.p)
ggplot() +
geom_boxplot(data = NA.month, aes(x = Month, y = Tave), alpha = 0.2, col = "darkblue", notch = TRUE) +
geom_boxplot(data = inst.month.T, aes(x = month, y = mean_temp_year), alpha = 0.2, col = "darkgreen", notch = TRUE) +
theme_bw() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
ggplot() +
geom_boxplot(data = NA.month, aes(x = Month, y = PPT), alpha = 0.2, col = "darkblue", notch = TRUE) +
geom_boxplot(data = inst.month.P, aes(x = month, y = total_precip_year), alpha = 0.2, col = "darkgreen", notch = TRUE) +
theme_bw() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
ggplot(data = months.T, aes(x=Year, y = Tave, col = Source)) +
geom_line() +
facet_grid(Month~.) +
theme_classic()
ggplot(data = months.P, aes(x=Year, y = PPT, col = Source)) +
geom_line() +
facet_grid(Month~.) +
theme_classic()
# ggplot() +
# geom_line(data = study_p, aes(x = month, y = total_precip), col= "skyblue", size = 1) +
# geom_ribbon(data = study_p, aes(x = month, ymin = (lower), ymax = (upper)), fill = "skyblue", alpha = 0.2) +
# geom_line(data = NA_p, aes(x = Month, y = mean_precip), col= "steelblue4", size = 1) +
# geom_ribbon(data = NA_p, aes(x = Month, ymin = (lower), ymax = (upper)), alpha = 0.2, fill= "steelblue4") +
# scale_x_continuous(breaks = c(1,2,3,4,5,6,7,8,9,10,11,12)) +
# theme_bw() +
# theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
# ggplot() +
# # London Data
# geom_line(data=london_temp_annual, aes(x = year, y = mean_temp), colour = "lightseagreen", linetype = "dashed", alpha = 0.8) +
# # Pinery Data
# geom_line(data=pinery_temp_annual, aes(x = year, y = mean_temp), colour="lightsalmon4", linetype = "dashed", alpha = 0.8) +
# # ClimateNA
# geom_line(data=annual, aes(x = Year, y = MAT), colour="darkgreen") +
# #Aesthetics
# xlab("Year") +
# ylab("Mean Annual Temperature (°C)") +
# scale_x_continuous(limits = c(1880, 2020), breaks=seq(1880,2020,20)) +
# #ggtitle("Average Annual Temperature in \n Southwestern Ontario from 1883 to 2017") +
# theme_bw() + theme(plot.title = element_text(hjust = 0.5)) +
# theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
i.m.p <- i.m.p %>% rename(Instrumental=PPT) %>% select(-Source)
m.m.p <- m.m.p %>% rename(ClimateNA=PPT) %>% select(-Source)
p.complete <- right_join(i.m.p, m.m.p, by = c("Year", "Month")) %>% na.omit()
p.matrix <- p.complete %>% select(-Year, -Month)
i.m.t <- i.m.t %>% rename(Instrumental=Tave) %>% select(-Source)
m.m.t <- m.m.t %>% rename(ClimateNA=Tave) %>% select(-Source)
t.complete <- right_join(i.m.t, m.m.t, by = c("Year", "Month")) %>% na.omit()
t.matrix <- t.complete %>% select(-Year, -Month)
m.t1 <- lm(data = t.complete, ClimateNA ~ Instrumental)
summary(m.t1)
##
## Call:
## lm(formula = ClimateNA ~ Instrumental, data = t.complete)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.41302 -0.35540 0.00878 0.33167 1.85683
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.418485 0.018835 22.22 <2e-16 ***
## Instrumental 0.978125 0.001511 647.44 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5205 on 1246 degrees of freedom
## Multiple R-squared: 0.997, Adjusted R-squared: 0.997
## F-statistic: 4.192e+05 on 1 and 1246 DF, p-value: < 2.2e-16
Remove outliers
p.complete <- p.complete %>% filter(Instrumental < 200)
Relevel
p.complete$Month <- relevel(p.complete$Month, ref = 5)
m.p1 <- lm(data = p.complete, ClimateNA ~ Instrumental)
summary(m.p1)
##
## Call:
## lm(formula = ClimateNA ~ Instrumental, data = p.complete)
##
## Residuals:
## Min 1Q Median 3Q Max
## -54.940 -11.719 -0.927 10.306 76.853
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 27.1542 1.2108 22.43 <2e-16 ***
## Instrumental 0.5411 0.0141 38.36 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 17.46 on 1234 degrees of freedom
## Multiple R-squared: 0.5439, Adjusted R-squared: 0.5435
## F-statistic: 1472 on 1 and 1234 DF, p-value: < 2.2e-16
m.p2 <- lm(data = p.complete, ClimateNA ~ Instrumental*Month)
summary(m.p2)
##
## Call:
## lm(formula = ClimateNA ~ Instrumental * Month, data = p.complete)
##
## Residuals:
## Min 1Q Median 3Q Max
## -47.608 -11.389 -1.025 10.188 79.701
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 19.680204 3.886351 5.064 4.74e-07 ***
## Instrumental 0.655126 0.043422 15.087 < 2e-16 ***
## Month1 10.222204 5.888583 1.736 0.082829 .
## Month2 2.832610 5.844245 0.485 0.627988
## Month3 2.557246 5.862190 0.436 0.662750
## Month4 5.253099 6.141415 0.855 0.392523
## Month6 11.500457 5.634142 2.041 0.041446 *
## Month7 21.570780 5.249084 4.109 4.23e-05 ***
## Month8 11.509542 5.343394 2.154 0.031439 *
## Month9 1.337119 5.507942 0.243 0.808231
## Month10 -1.132427 5.649630 -0.200 0.841168
## Month11 -0.894239 6.150210 -0.145 0.884419
## Month12 15.112783 6.166552 2.451 0.014396 *
## Instrumental:Month1 -0.180905 0.066507 -2.720 0.006620 **
## Instrumental:Month2 -0.096321 0.076943 -1.252 0.210869
## Instrumental:Month3 -0.062959 0.076407 -0.824 0.410105
## Instrumental:Month4 0.006738 0.073945 0.091 0.927413
## Instrumental:Month6 -0.081111 0.065660 -1.235 0.216955
## Instrumental:Month7 -0.262228 0.058347 -4.494 7.65e-06 ***
## Instrumental:Month8 -0.204063 0.059557 -3.426 0.000632 ***
## Instrumental:Month9 -0.029312 0.060796 -0.482 0.629792
## Instrumental:Month10 -0.017195 0.064328 -0.267 0.789279
## Instrumental:Month11 -0.021626 0.066670 -0.324 0.745708
## Instrumental:Month12 -0.266745 0.066686 -4.000 6.72e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.84 on 1212 degrees of freedom
## Multiple R-squared: 0.5832, Adjusted R-squared: 0.5753
## F-statistic: 73.73 on 23 and 1212 DF, p-value: < 2.2e-16
t.complete$fit.1 <- fitted(m.t1) #Model predictions or fitted values
p.complete$fit.1 <- fitted(m.p1)
p.complete$fit.2 <- fitted(m.p2)
ggplot(data = t.complete) +
geom_line(aes(x=Instrumental, y = fit.1, col = Month)) +
geom_point(aes(x=Instrumental, y = ClimateNA, col = Month), size = 0.5, alpha = 0.3) +
labs(x = "Instrumental: Average Monthly Temperature",
y = "ClimateNA: Average Monthly Temperature") +
theme_classic()
Relevel to plot
p.complete$Month <- factor(p.complete$Month, levels = c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "11", "12"))
ggplot(data = p.complete) +
geom_line(aes(x=Instrumental, y = fit.1, col = Month)) +
geom_point(aes(x=Instrumental, y = ClimateNA, col = Month), size = 0.5, alpha = 0.5) +
labs(x = "Instrumental: Monthly Precipitation (mm)",
y = "ClimateNA: Monthly Precipitation (mm)") +
theme_classic()
ggplot(data = p.complete) +
geom_line(aes(x=Instrumental, y = fit.2, col = Month)) +
geom_point(aes(x=Instrumental, y = ClimateNA, col = Month), size = 0.5, alpha = 0.4) +
lims(x = c(0,200), y = c(0,200)) +
labs(x = "Instrumental: Monthly Precipitation (mm)",
y = "ClimateNA: Monthly Precipitation (mm)") +
theme_classic()
Normality:
plot(m.t1, which=2)
plot(m.p1, which=2)
plot(m.p2, which=2)
Equal variance:
plot(m.t1, which=3)
plot(m.p1, which=3)
plot(m.p2, which=3)
Linearity:
plot(m.t1, which=1)
plot(m.p1, which=1)
plot(m.p2, which=1)
Leverage:
plot(m.t1, which=4)
plot(m.p1, which=4)
plot(m.p2, which=4)
Sys.time()
## [1] "2020-06-22 16:32:42 PDT"
git2r::repository()
## Local: master /Users/JenBaron/Documents/UWO/UWO NSERC/Statistical Analysis/Climate/Climate_Pinery
## Remote: master @ origin (https://github.com/JenBaron/Climate_Pinery.git)
## Head: [ad4302f] 2020-06-19: seasonal climateNA
sessionInfo()
## R version 4.0.0 (2020-04-24)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Catalina 10.15.5
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_CA.UTF-8/en_CA.UTF-8/en_CA.UTF-8/C/en_CA.UTF-8/en_CA.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] Hmisc_4.4-0 Formula_1.2-3 survival_3.1-12 lattice_0.20-41
## [5] gridExtra_2.3 tidyr_1.0.2 dplyr_1.0.0 ggplot2_3.3.1
##
## loaded via a namespace (and not attached):
## [1] tidyselect_1.1.0 xfun_0.14 purrr_0.3.4
## [4] splines_4.0.0 colorspace_1.4-1 vctrs_0.3.1
## [7] generics_0.0.2 htmltools_0.4.0 yaml_2.2.1
## [10] base64enc_0.1-3 rlang_0.4.6 pillar_1.4.4
## [13] foreign_0.8-79 glue_1.4.1 withr_2.2.0
## [16] RColorBrewer_1.1-2 jpeg_0.1-8.1 lifecycle_0.2.0
## [19] stringr_1.4.0 munsell_0.5.0 gtable_0.3.0
## [22] htmlwidgets_1.5.1 evaluate_0.14 labeling_0.3
## [25] latticeExtra_0.6-29 knitr_1.28 htmlTable_1.13.3
## [28] Rcpp_1.0.4.6 acepack_1.4.1 scales_1.1.1
## [31] backports_1.1.7 checkmate_2.0.0 farver_2.0.3
## [34] png_0.1-7 digest_0.6.25 stringi_1.4.6
## [37] grid_4.0.0 tools_4.0.0 magrittr_1.5
## [40] tibble_3.0.1 cluster_2.1.0 crayon_1.3.4
## [43] pkgconfig_2.0.3 ellipsis_0.3.1 Matrix_1.2-18
## [46] data.table_1.12.8 rmarkdown_2.1 rstudioapi_0.11
## [49] R6_2.4.1 rpart_4.1-15 git2r_0.26.1
## [52] nnet_7.3-14 compiler_4.0.0